Low- Temperature Phase Boundary of dilute Lattice Spin Glasses 
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The thermal-to-percolative crossover exponent 0, well-known for ferromagnetic systems, is studied 
extensively for Ed wards- Anderson spin glasses. The scaling of defect energies are determined at the 
bond percolation threshold pc, using an new algorithm. Simulations extend to system sizes above 
A'^ = 10* in dimensions d = 2, . . . , 7. The results can be related to the behavior of the transition 
temperature Tg ^ {p ~ pc)"^ between the paramagnetic and the glassy regime for p \ p^. In three 
dimensions, where our simulations predict (j) — 1.127(5), this scaling form for Tg provides a rare 
experimental test of predictions arising from the equilibrium theory of low-temperature spin glasses. 
For dimension near and above the upper critical dimension, the results provide a new challenge to 
reconcile mean-field theory with finite-dimensional properties. 

PACS numbers: 75.10.Nr , 02.60.Pn , 64.70.Pf 



The exploration of low-temperature properties of dis- 
ordered systems remains an important and challenging 
problemjii^ The paradigmatic model for such phenom- 
ena is the Edwards- Anderson (EA) spin glass, 



(1) 



Disorder effects arise via quenched random bonds, J^j , 
mixing fcrro- and anti-ferromagnetic couplings between 
nearest-neighbor spins, that lead to conflicting con- 
straints and frustrated variables. It is believed that an 
understanding of static and dynamic features of EA may 
aid a description of the unifying principles expressed in 
glassy materials.^ Most insights into finite-dimensional 
systems has been gained through computational ap- 
proaches that elucidate low-T properties ^'^'-S 

Here, we extract the response induced through defect- 
interfacesii^ at T = 0, created by fixing the spins along 
the two faces of the open boundary in one direction. 
Ground state energies Eq and Eq of an instance of size 
N = L'^ are determined that differ by reversing all spins 
on one of the faces. The distribution P{AE) of interface 
energies AE = Eq — Eo created by this perturbation of 
scale L on the boundary is obtained. The typical energy 
scale, represented by the deviation a{AE), grows as 



a{AE) - 



(2) 



This relation defines^iii^ the stiffness exponent y charac- 
terizing the defect energy, a fundamental quantity assess- 
ing low-temperature fluctuations: a positive value of y, 
as found^ in EA for d > di ~ 5/2, denotes the increase 
in the energetic cost accompanying a growing number 
of variables perturbed from their position in the ground 
state (i. e., "stiffness"). The rise in strain for stronger 
disturbances signals the presence of an ordered state. In 
turn, for systems with y < such order is destabilized 
by fluctuations that spread unimpeded. 

Instead of determining the interface scaling on a com- 
pact lattice structure, we will focus here on the interface 
energy a{AE) on a bond-diluted lattice, in particular, at 
the percolation threshold Pc, see Fig. [TJ Due to the ten- 
uous fractal nature of the percolating cluster at Pc, no 



long-range order can be sustained, andi^ 

<j{AE)l,p^ ~ L^^ with yp < 0, 



(3) 



i. e., defects possess a vanishing interface energy. In- 
terest in the exponent yp stems from its relation to 
the "thermal-percolative cross-over exponent" (p defined 



Tgip) [p-PcY 



with 



-vyp 



(4) 



where v is the correlation-length exponent associated 
with lattice percolation, ^^'^^ ^ ^ {p — Pc)~'^ ■ Of par- 
ticular experimental interest is the result for d = 3, 
yp = -1.289(6), predicting ^ 1.127(5) with ly = 
0.87436(46)." All resuhs for d = 2, . . . , 7 are listed in 
Tab.U 

The exponent has been studied intensely numerically, 
theoretically, and experimentall y^^'^^i^'^i^^i^^i^°i^^i^^ for 
ferromagnetic systems some 30 years ago, and just re- 
cently was discussed for quantum spins. But aside from 
its initial treatment in Ref.[l3, there are no investigations 
on spin glasses. This is even more surprising, since this 
exponent provides a non-trivial, experimentally testable 
prediction derived from scaling arguments of the equi- 
librium theory at low temperatures. Such tests are few 
as disordered materials by their very nature fall out of 
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FIG. 1: (Color online) Phase diagram for bond-diluted spin 
glasses {d > di). In the spin-glass phase (SG) for T < Tg and 
p > Pc y in Eq. ((2]) is > 0, while y ^ yp < in Eq. Q at 
p = Pc and T — Tg — 0. In the paramagnetic phase (PM) for 
p < Pc defects decay exponentially for all T. The exponent 
(f) in Eq. Q describes the boundary Tg{p) for p\pc- 
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TABLE L List of the parameters used and exponents found 
in our simulations for d = 2, . . . ,7 . Lmax denotes the largest 
lattice size considered. We have used the bo nd-p ercolation 
thresholds Pc from Ref. [3 for d = 3 and Ref. [iSl for d > 4. 
The correlation-length exponents u for percolation are from 
Ref. [l^ in d = 3 and from Ref. [ll| for d > 4, where = 1/2 is 
exact above the upper critical dimension, d > 6. 
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equilibrium when entering the glassy state. The phase 
boundary itself provides the perfect object for such a 
study: It can be approached by theory from below and by 
experiments from above where equilibration is possible. 

There is reason to believe that the phase bound- 
ary in Eq. ^ and Fig. [T] is experimentally accessi- 
ble for certain materials. Ref. |2J already provided 
highly accurate results for the freezing temperature Tm 
as a function of dilution x for a doped, cristalline 
glass, (Lai_a;Gda;)8oAu2o, proposing a linear dependence, 
Tm ~ x. The tabulated data is equally well fitted by 
Eq. ([U in that regime. Ref. [H determined a phase dia- 
gram for (Fea:Nii_a:)75Pi6B6Al4, an amorphous alloy, for 
a wide range of temperatures T and site-concentrations 
X but did not discuss its near-linear behavior at low x. A 
similar phase diagram for the insulator CdCr2£!;In2(i_a;)S4 
can be found in Fig. 1.1a of Ref. [2^. New experiments 
dedicated to the limit x \, Xc should provide results of 
sufficient accuracy to test our prediction for 

A match of computational prediction and experiment 
would lend credibility to the EA model and its simpli- 
fying assumptions, such as universality with respect to 
the details of the bond distribution P{J), an issue re- 
cently revisited by Ref. [27l. Our simulations, conducted 
here for Gaussian bonds, can be repeated for any P{J) of 
zero mean and unit variance, but would significantly in- 
crease computational cost. Simply to demonstrate that 
theoretically such universality exists, we have repeated 
our simulations with a Lorcntzian bond distribution ar- 
gued for by Ref. [28l, and with powcrlaw-distributed bonds 
P{J) cx |J|"-i for |J| > 1 at a -1/2. This compari- 
son, presented below in Fig. [31 clearly show reproducibil- 
ity for a wide class of P{J)- Ref. Tc? also considered 
power-law distributions, but for a > and |J| < 1, to 
point out that the response to perturbations is in prin- 
ciple non- universal at Pc'- Without long-range order, it 
is y < in Eq. ^ and energy scales are not diverg- 
ing. Then, y (and cj)) become dependent on the details 
of P{J) near J = 0, and Ref. 10 finds an interesting 
change in behavior for a < ac ~ 0.75. Such a diverg- 
ing bond distribution results from integrating the RKKY 



couplings over many weak bonds in its far-distance tail. 
In realistic materials, such bonds are screened out (see, 
for instance, Ref. [13), and P{J) is bounded, justifying 
the use of Gaussian bonds, say. 

Following the discussion in Refs. lToll29l for diluted lat- 
tices at p — > we have to generalize the scaling relation 
for the defect energy (t{AE) in Eq. ^ to 

a{AE)L,p y{p)Lyf {L/i{p)) . (5) 

Here, y ^ [p ~ PcY ~ ^^*-/'^ is surface tension and 
~ (p ~ Pc)~'^ is the correlation length for perco- 
lation. The scaling function / is defined to be constant 
for L ^ £^{p) ^ 1, where percolation (and hence, ^) plays 
no role and we regain Eq. ^ for p > pc- 

For ^ 3> i 1, Eq. ([5]) requires f{x) ~ x'^ for x — > 
to satisfy ct ^ with some power of L, needed to cancel 
the ^-dependence 'At p — pc- Thus, /i = —t/v^ and if we 
define yp = y + p = y — t/v to mark the L-dependence of 
cr at p = Pc as in Eq. ([3]), we get t = v[y — yp). Finally, at 
the cross-over ^ ~ L, where the range L of the excitations 
cr{AE) reaches the percolation length beyond which spin 
glass order ensues, Eq. ^ yields with from Eq. ([3]), 

a(Ai?)^(,)^, ~ (p - Pc)* - (P ~ Pet ■ (6) 

Associating a temperature with this cross-over by 
(t(A_E)j(p) p ~ Tg (for T > Tg, thermal fluctuation de- 
stroy order), leads to Eq. ([4]), relating p and Tg. 

In our simulations we have used the method of bond re- 
ductions described previouslyi^i^'^ A set of rules is de- 
fined and applied recursively to trace out spins assuming 
that T = 0. These exact rules apply to general Ising spin 
glass Hamiltonians as in Eq. ([T]) with any bond distri- 
bution P{J), discrete or continuous, on arbitrary sparse 
graphs, and lead to fewer but more highly interconnected 
spins and renormalized bonds, see Ref. 1321. Starting from 
a Hamiltonian as in Eq. ([1]), in general, new terms are 
generated by this procedure that have not been part of 
the Hamiltonian before, such as multi-spin interaction. 
Although the number of spins decreases one-by-one, the 
number of new terms grows exponentially and the pro- 
cedure usually becomes inefficient. Yet, near pc, we can 
apply a subset of these rules efficiently while leaving the 
form of the 2-spin Hamiltonian in Eq. ([T]) invariant. 

Our recursive set of rules is based on the following ob- 
servations. Near pc, most spins have a low degree of inter- 
connectivity; on average, that degree fluctuates around 
unity in any dimension d. In fact, many spins are entirely 
disconnected, do not contribute to the Hamiltonian, and 
thus can be discarded. Degree-1 spins can always be sat- 
isfied and are easily traced out, with their bond weight 
always (at T = 0) lowering the energy. Once all degree-1 
spins have been recursively traced out, any degree-2 spin 
can be reduced also by replacing it by a new bond be- 
tween its two neighbors and another off-set to the global 
energy. Having reduced all degree-1 and -2 spins, there is 
even a "star-triangle" rule to reduce any degree-3 spins 
while only producing new 2-spin interactions between its 
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neighbors. Although this step in principle could create 
a 3-spin interaction not present in the Hamiltonian in 
Eq. ^ , all such terms involving an odd number of spins 
vanish due to the Z2-symmetry of the Ising spins. 

A new rule^ that proved particularly effective at pc fo- 
cuses on spins of arbitrary degree but with "superbonds" . 
A spin Xi has a superbond, if one bond's absolute weight 
dominates, | Ji,fc| > Sj^^fc l^^ij li other bonds attached 
to Xi. In the ground state (T = 0), that bond is always 
satisfied and its spin determined by its neighbor along 
that bond. This rule often triggers new avalanches of 
further reductions with the simpler rules. 

Previously, we have applied these rules above pc to 
study the defect energy within the spin glass state (SG), 
see Fig. [T] Considering dilute lattices with p > Pc but 
well below p = 1 allowed the study of larger lattice sizes 
L for improved scaling, and produced results^i^ in di- 
mensions up to d = 7, unattainable with undiluted lat- 
tices. For p > Pc an optimization heuristic was essential 
to approximate the ground state of the remainder graph, 
consisting of highly interconnected spins that remain af- 
ter all reduction rules have been exhausted. In contrast, 
at p = Pc, these remained graphs are - almost - gone 
entirely. Thus, the attainable system sizes N = L'^ are 
nearly unrestricted and have reached well above N — 10^ 
in our simulations, mostly limited by the need to gen- 
erate sufficient statistics (i. e., about 10'* instances for 
N = 25^ or 15'^). Yet, in d = 2 and 3, the remain- 
der graphs are the limiting factor on system sizes (at 
about N = 10^). Although remainders have less than 100 
spins, typically well-approximated with a good heuris- 
tic, we implemented costly exact methods^ to optimize 
them. The slightest inaccuracy affected the statistical 
averages, as defect energies AE are the difference of two 
almost equal ground state energies Eq and Eq. One tech- 



nical problem in implementing our algorithm with such 
large system sizes is posed by memory limitations. In- 
stead of constructing an entire lattice with N spins, each 
with potentially 2d bonds, before applying the reduction 
rules, we build up the L'^-spin lattice as a sequence of 
L hyper-planes of L"^^^ spins. During the process, we 
keep the first and the most recently added plane fixed, 
but already reduce recursively all spins in the interven- 
ing planes as far as possible, before the next hyper-plane 
is added. This process requires extensive bookkeeping 
and backtracking which can be done fast while reducing 
memory use by ~ 1/L. 

In Fig. [2] we present all data of our simulations for 
d = 2, . . . , 7 in an extrapolation plot. In Fig. [3] we com- 
pare the same data from d = 3 (d = 4) together with 
those from the Lorentzian (power-law) bond distribution 
P(J), as discussed above. Since our data reaches above 
the upper critical dimension = 6 (of both, percola- 
tion and spin glasses) and should approach mean-field 
behavior, it is most natural to replace L^^ with N^/'^ 
in Eq. ([3]) and extrapolate for yp/d ~ ln((T)/ In(A^). As 
Fig. [2] shows, aside from d = 2 < di, the extrapolations 
for d > di all seem to share common characteristics and 
appear to vary smoothly with d. In particular, we have 
pushed the simulations in d = 7 to large enough N to con- 
clude that there appears to be no drastic change in the 
scaling behavior above du. Increasing corrections make it 
harder to reach asymptotic scaling beyond d = 7. Inter- 
estingly, all transients in Fig. O themselves extrapolate to 
an intercept consistent with —1/2, indicative of a higher- 
order correction term with a d- independent exponent. 

The result in d = 2, where Tg = 0, is very close to 
that theoretically predicted in Ref. [13, yp « —0.99, and 
could conceivably be = — 1 exactly, as it is in d = 1 
(where pc = 1). That would suggest that the spin glass 
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FIG. 2: (Color Online) Plot of a (AE) in Eq. © as a func- 
tion of system size N = L'' in extrapolated form. Plotting 
ln{a)/ln{N) vs. l/ln{N) and linearly extrapolating (dashed 
lines), we extract the asymptotic values for yp/d in Tab.|T]at 
1 / ln(A'^) = 0. Note that the fitted asymptotic regime here 
corresponds to orders of magnitude of scaling in a (less in- 
sightful) plot of Incr vs. In A'^. Note the increasing corrections 
to scaling for larger d before asymptotic behavior is obtained. 
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FIG. 3: (Color Online) Comparison of the data for d = 3 (top) 
and d = 4 (bottom) for different bond distributions P{J), 
plotted in the same way as Fig. (2] There is little difference 
between the Gaussian and the Lorentzian bonds even in d = 3, 
as both are similarly smooth near J = 0. For our power-law 
bonds with vanishing support for jjj < 1, our methods are 
very inefficient. Yet, at least in d = 4 those bonds produce 
results far from but consistent with Gaussian bonds. 
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on a d = 2 percolation cluster essentially consists of an 
extensive linear backbone of bonds. 

In Fig. m we plotted yp/d and (f) from Tab. H] vs. 
1/d to explore the large-d limit. This extrapolation plot 
suggests a trend towards a vanishing value for yp/d at 
l/d — 0, i. e., yp varies sub-linearly with d. In compar- 
ison, the data for the stiffness exponents y/d inside the 
spin-glass regime replotted from Ref. 9 appears consis- 
tent with the prediction^i^ of y/d ~ 1/6. It would be 
difficult to suspect a systematic bias in the apparent drift 
of the high-d data points for yp/d, as the computations 
are exact. Yet, statistical errors clearly become increas- 
ingly significant for larger d, see Fig. [2] It is not obvious 
how to directly obtain yp for d = oo, which may corre- 
spond to a (replica-symmetric) T = Tg = Viana-Bray 
modeP^ at the Erdos-Renyi percolation point. [Such a 
calculation has been undertaken for the fully connected 
(replica-symmetry broken, T < Tg) SK modeli^l] Finally, 
we note a distinct minimum in cf), with (/)6 ~ 1, exactly at 
the upper critical dimension d„ = 6, due to the product 
of increasing \yp\ and decreasing v. 
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FIG. 4: (Color online) Plot of the exponents yp/d (top) and 
phi (bottom) in Tab. [T] as a function of l/d. For comparison, 
also plotted (top) are the stiffness exponents y/d inside the 
spin-glass regime^^ (p > pc). Some of the large error bars for 
(j) = —vyp originate with uncertainties in v^. 
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